Zanim przejdziemy do wyprowadzenie czterech fundamentalnych równań struktury gwiazd, omówimy prosty model analityczny. Cztery docelowe równania, w uproszczonej formie, która wystarczana potrzeby tego wykładu, to:

1. Równanie ciągłości D[m[r],t]==4 Pi r^2 ρ[r]
2. Równanie równowagi hydrostatycznej: D[P[r],r] == -G m[r] ρ[r]/r^2
3. Równanie transportu energii (przez promieniowanie lub konwenkcję): D[T[r],r] == ?
4. Równanie produkcji energii: D[L[r],r] == ?

Cztery niewiadome funkcje: m[r] (rozkład masy), P[r] (rozkład ciśnienia) , T[r] (rozkład temperatury) , L[r] (jasność gwiazdy) można w sensowny sposób przybliżyć analitycznie. Jednym z takich ugruntowanych w astrofizyce modeli jest m. Eddingtona. 

Jego założenia to:

a) stosunek ciśnienia promieniowania i gazu doskonałego wszędzie jest taki sam => równanie stanu okazuje się politropowe z n=3. 
b) cała produkowana energia pojawia się w cudowny sposób znikąd w r=0, co daje L[r] = const = (Subscript[L, ⊙])

W modelu Eddingtona wszystkie wielkości wyrażają się przez funkcję Lane-Emdena (Subscript[θ, 3]).

(*Solve[{Pgaz==kB T/mu/μ ρ,Prad==1/3 a T^4,Pgaz/P==β,Prad+Pgaz==P},{P,Prad,Pgaz,T}]*)

Idea modelu Eddingtona to miks gazu fotonowego i doskonałego w stałej proporcji, tak aby stosunek ciśnień 2 składników był stały

 .

$Assumptions={kB>0,μ>0,0<β<1,mu>0}

Dwa składniki to: 
	1) gaz fotonowy (promieniowanie) o równaniu stanu ; ( - gęstość energii gazu fotonowego,  - stała promieniowania,  - temperatura)
	2) gaz doskonały ;     ( - stała Boltzmana,  - średnia waga molekularna cząstek gazu, np:  dla wodoru neutralnego,  - atomowa jednostka masy)

Eliminate[{Pgaz==kB T/mu/μ ρ,Prad==1/3 a T^4,Pgaz/P==β,Prad+Pgaz==P},{Prad,Pgaz,T}]

%/.P^3->P3

Solve[%,P3]

K=(P3/.%[[1]]//FullSimplify)^(1/3)/.ρ->1

Okazuje się, że , czyli równanie stanu jest politropowe z wykładnikiem , czyli indeksem politropowym . Masa jest więc masą Chandrasekhara.

M=Sqrt[16 K^3/Pi/G^3](-x0^2 x1)

(*Solve[M==Msun,β]*)

kBDef=Quantity[1,"BoltzmannConstant"]//UnitConvert//QuantityMagnitude

kB=N[kBDef]

c=299792458;

mu=Quantity[1,"AtomicMassUnit"]//UnitConvert//QuantityMagnitude//N

aDef=Quantity[1,"RadiationConstant"]//UnitConvert//QuantityMagnitude

a=N[aDef]

Import["https://th.if.uj.edu.pl/~odrzywolek/homepage/students/aa/Lane-Emden.m"]

x0=zero[3]

x1=slope[3]

(*x0=6.8968486193769614492810885206340348229744082117093697764206*)

(*x1=-0.0424297576044513703852644025245080863173417829990589903735*)

μ=16/27

G=Quantity[1,"GravitationalConstant"]//UnitConvert//QuantityMagnitude//N

Msun=Quantity[1,"SolarMass"]//UnitConvert//QuantityMagnitude//N

βSun=Solve[M==Msun,β,Reals][[1]]

1-β/.%

M/.βSun

18.1/μ^2

Model Eddingotona sensownie przewiduje spodziewany zakres mas gwiazd, od ok. 0.01  do 100  .

LogLogPlot[M/Msun/.β->1-beta,{beta,10^-8,1},PlotRange->{0.001,1000},AxesLabel->{"1-β=(Subscript[P, rad])/P","M/(Subscript[M, ⊙])"},GridLines->{{1},{1}}]

LogLogPlot[M/Msun/.β->beta,{beta,10^-1,1},PlotRange->{0.001,1000},AxesLabel->{"β=(Subscript[P, gaz])/P","M/(Subscript[M, ⊙])"},GridLines->{{1},{1}}]

M/Msun/.β->0.5 (* Masa gwiazdy w której ciśnienia promieniowania i materii zrównują się *)

R=Quantity[1,"SolarRadius"]//UnitConvert//QuantityMagnitude//N

Lsun=Quantity[1,"SolarLuminosity"]//UnitConvert//QuantityMagnitude//N

Tsun=Quantity[1,"SolarEffectiveTemperature"]//UnitConvert//QuantityMagnitude